/**
   simulator.cu - Created by Timothy Morey on 9/2/2011.
 */

#include <stdio.h>
#include <string>
using namespace std;

#include <boost/timer.hpp>

#include "device_controller.h"
#include "glut_display.h"
#include "parts.cuh"
#include "renderer.h"
#include "simulator.h"
#include "../common/vis_server.h"

Simulator* Simulator::_Instance = 0;

Simulator::Simulator()
  : _HostX(0), _HostY(0), 
    _HostVX(0), _HostVY(0), 
    _HostAX(0), _HostAY(0), 
    _HostFrame(0),
    _FrameWidth(512),
    _FrameHeight(512),
    _PauseFlag(false),
    _StopFlag(false),
    _Vis(0)
{
  _Parameters.N = 1000;
  _Parameters.MinX = 0.0;
  _Parameters.MaxX = 100.0;
  _Parameters.MinY = 0.0;
  _Parameters.MaxY = 100.0;
  _Parameters.dt = 0.001;
  
  size_t memsize = _Parameters.N * sizeof(DoF);
  CUDA_CALL( cudaHostAlloc(&_HostX, memsize, cudaHostAllocDefault) );
  CUDA_CALL( cudaHostAlloc(&_HostY, memsize, cudaHostAllocDefault) );
  CUDA_CALL( cudaHostAlloc(&_HostVX, memsize, cudaHostAllocDefault) );
  CUDA_CALL( cudaHostAlloc(&_HostVY, memsize, cudaHostAllocDefault) );
  CUDA_CALL( cudaHostAlloc(&_HostVM, memsize, cudaHostAllocDefault) );
  CUDA_CALL( cudaHostAlloc(&_HostAX, memsize, cudaHostAllocDefault) );
  CUDA_CALL( cudaHostAlloc(&_HostAY, memsize, cudaHostAllocDefault) );

  size_t framesize = _FrameWidth * _FrameHeight * 4 * sizeof(unsigned char);
  CUDA_CALL( cudaHostAlloc(&_HostFrame, framesize, cudaHostAllocDefault) );
  
  memset(_HostX, 0, memsize);
  memset(_HostY, 0, memsize);
  memset(_HostVX, 0, memsize);
  memset(_HostVY, 0, memsize);
  memset(_HostVM, 0, memsize);
  memset(_HostAX, 0, memsize);
  memset(_HostAY, 0, memsize);
  memset(_HostFrame, 0, framesize);
  
  this->RandomizePositions();
}

Simulator::~Simulator()
{
  // TODO: stop devices if they are running

  CUDA_CALL( cudaFreeHost(_HostX) );
  CUDA_CALL( cudaFreeHost(_HostY) );
  CUDA_CALL( cudaFreeHost(_HostVX) );
  CUDA_CALL( cudaFreeHost(_HostVY) );
  CUDA_CALL( cudaFreeHost(_HostVM) );
  CUDA_CALL( cudaFreeHost(_HostAX) );
  CUDA_CALL( cudaFreeHost(_HostAY) );
  CUDA_CALL( cudaFreeHost(_HostFrame) );
}

Simulator& Simulator::Instance()
{
  if(! _Instance)
    _Instance = new Simulator();

  return *_Instance;
}

void Simulator::DestroyInstance()
{
  if(_Instance)
  {
    delete _Instance;
    _Instance = 0;
  }
}

bool Simulator::UpdateParameters(const Parameters& parms)
{
  bool success = true;

  bool resume = false;
  if(this->IsRunning() && ! _PauseFlag)
  {
    this->Pause();
    resume = true;
  }

  if(parms.N != _Parameters.N)
  {
    size_t memsize = parms.N * sizeof(DoF);

    DoF* newHostX = 0;
    DoF* newHostY = 0;
    DoF* newHostVX = 0;
    DoF* newHostVY = 0;
    DoF* newHostVM = 0;
    DoF* newHostAY = 0;
    DoF* newHostAX = 0;
    DoF* oldHostX = _HostX;
    DoF* oldHostY = _HostY;
    DoF* oldHostVX = _HostVX;
    DoF* oldHostVY = _HostVY;
    DoF* oldHostVM = _HostVM;
    DoF* oldHostAX = _HostAX;
    DoF* oldHostAY = _HostAY;

    CUDA_CALL( cudaHostAlloc(&newHostX, memsize, cudaHostAllocDefault) );
    CUDA_CALL( cudaHostAlloc(&newHostY, memsize, cudaHostAllocDefault) );
    CUDA_CALL( cudaHostAlloc(&newHostVX, memsize, cudaHostAllocDefault) );
    CUDA_CALL( cudaHostAlloc(&newHostVY, memsize, cudaHostAllocDefault) );
    CUDA_CALL( cudaHostAlloc(&newHostVM, memsize, cudaHostAllocDefault) ;)
    CUDA_CALL( cudaHostAlloc(&newHostAX, memsize, cudaHostAllocDefault) );
    CUDA_CALL( cudaHostAlloc(&newHostAY, memsize, cudaHostAllocDefault) );

    memset(newHostX, 0, memsize);
    memset(newHostY, 0, memsize);
    memset(newHostVX, 0, memsize);
    memset(newHostVY, 0, memsize);
    memset(newHostVM, 0, memsize);
    memset(newHostAX, 0, memsize);
    memset(newHostAY, 0, memsize);

    for(int i = 0; i < _Parameters.N && i < parms.N; i++)
    {
      newHostX[i] = _HostX[i];
      newHostY[i] = _HostY[i];
      newHostVX[i] = _HostVX[i];
      newHostVY[i] = _HostVY[i];
      newHostVM[i] = _HostVM[i];
      newHostAX[i] = _HostAX[i];
      newHostAY[i] = _HostAY[i];
    }

    int oldN = _Parameters.N;

    if(oldN > parms.N)
      _Parameters.N = parms.N;
    _HostX = newHostX;
    _HostY = newHostY;
    _HostVX = newHostVX;
    _HostVY = newHostVY;
    _HostVM = newHostVM;
    _HostAX = newHostAX;
    _HostAY = newHostAY;
    if(oldN < parms.N)
      _Parameters.N = parms.N;

    CUDA_CALL( cudaFreeHost(oldHostX) );
    CUDA_CALL( cudaFreeHost(oldHostY) );
    CUDA_CALL( cudaFreeHost(oldHostVX) );
    CUDA_CALL( cudaFreeHost(oldHostVY) );
    CUDA_CALL( cudaFreeHost(oldHostVM) );
    CUDA_CALL( cudaFreeHost(oldHostAX) );
    CUDA_CALL( cudaFreeHost(oldHostAY) );

    if(_Parameters.N > oldN)
      success = this->RandomizePositions(oldN, _Parameters.N - 1);
  }
  
  if(parms.MinX != _Parameters.MinX ||
     parms.MaxX != _Parameters.MaxX ||
     parms.MinY != _Parameters.MinY ||
     parms.MaxY != _Parameters.MaxY)
  {
    _Parameters.MinX = parms.MinX;
    _Parameters.MaxX = parms.MaxX;
    _Parameters.MinY = parms.MinY;
    _Parameters.MaxY = parms.MaxY;
    
    for(int i = 0; i < _Parameters.N; i++)
      EnsureInBounds(_HostX[i], _HostY[i], 
                     _Parameters.MinX, _Parameters.MaxX,
                     _Parameters.MinY, _Parameters.MaxY);

    GlutDisplay::Instance().UpdateScreenBounds();
  }

  if(resume)
    this->Resume();

  return success;
}

void Simulator::SetStepSize(DoF dt)
{
  _Parameters.dt = dt;
}

bool Simulator::RandomizePositions(int first /* = -1 */, int last /* = -1 */)
{
  bool retval = true;
  if(_HostX && _HostY)
  {
    srand(time(0));

    if(first == -1)
      first = 0;

    if(last == -1)
      last = _Parameters.N - 1;
    
    for(int i = first; retval && i <= last; i++)
    {
      bool acceptable = false;
      int numAttempts = 0;

      while(! acceptable)
      {
	_HostX[i] = _Parameters.MinX + 
	  (_Parameters.MaxX - _Parameters.MinX) * ((DoF)rand() / (DoF)RAND_MAX);
	_HostY[i] = _Parameters.MinY + 
	  (_Parameters.MaxY - _Parameters.MinY) * ((DoF)rand() / (DoF)RAND_MAX);
	
	acceptable = true;
	for(int j = 0; j < _Parameters.N; j++)
	{
	  if(j == i)
	  {
	    j = last;
	  }
	  else
	  {
	    if(DistanceTo(_HostX[i], _HostY[i], 
			  _HostX[j], _HostY[j],
                          _Parameters.MinX, _Parameters.MaxX,
                          _Parameters.MinY, _Parameters.MaxY) < 1.0)
	    {
	      acceptable = false;
	      break;
	    }
	  }
	}

	numAttempts++;

	if(numAttempts > 1024)
	{
	  // Then we're having trouble placing the particle, and it is likely
	  // that there really isn't enough room in the world for the requested
	  // number of particles.  Avoid infinite loop.

	  retval = false;
	  break;
	}
      }
    }
  }

  return retval;
}

void Simulator::RandomizeVelocities(DoF magnitude /* = 1.0 */,
				    int first /* = -1*/, 
				    int last /* = -1 */)
{
  if(_HostVX && _HostVY)
  {
    srand(time(0));

    if(first == -1)
      first = 0;

    if(last == -1)
      last = _Parameters.N - 1;

    for(int i = first; i <= last; i++)
    {
      _HostVX[i] = (-magnitude) + (2 * magnitude) * ((DoF)rand() / (DoF)RAND_MAX);
      _HostVY[i] = (-magnitude) + (2 * magnitude) * ((DoF)rand() / (DoF)RAND_MAX);
      _HostVM[i] = sqrt(_HostVX[i] * _HostVX[i] + _HostVY[i] * _HostVY[i]);
    }
  }
}

void Simulator::SetFrameSize(unsigned int width, unsigned int height)
{
  printf("SetFrameSize!\n");
  _FrameWidth = width;
  _FrameHeight = height;
  size_t framesize = _FrameWidth * _FrameHeight * 4 * sizeof(unsigned char);

  unsigned char* newFrame = 0;
  CUDA_CALL( cudaHostAlloc(&newFrame, framesize, cudaHostAllocDefault) );
  memset(newFrame, 0, framesize);

  // TODO: if rendering is in progress on another thread, deleting the host
  // memory may crash.  Need some way to warn users of the frame buffer.
  CUDA_CALL( cudaFreeHost(_HostFrame) );

  _HostFrame = newFrame;

  if(_Vis)
    _Vis->SetImage(_HostFrame, _FrameWidth, _FrameHeight);
}

DoF Simulator::GetX(int i) const
{
  DoF retval = 0.0f;

  if(_HostX && i >= 0 && i < _Parameters.N)
    retval = _HostX[i];

  return retval;
}

DoF Simulator::GetY(int i) const
{
  DoF retval = 0.0f;

  if(_HostY && i >= 0 && i < _Parameters.N)
    retval = _HostY[i];

  return retval;
}

DoF Simulator::GetVM(int i) const
{
  DoF retval = 0.0f;

  if(_HostVM && i >= 0 && i < _Parameters.N)
    retval = _HostVM[i];

  return retval;
}

void Simulator::PrintStepRateReport() const
{
  std::vector<DeviceController*>::const_iterator iter;
  printf("%6s | %14s | %14s | %8s | %16s\n", "Device", "CPU Inner Time", 
         "CPU Outer Time", "GPU Time", "Steps Per Second");
  printf("%6s | %14s | %14s | %8s | %16s\n", "------", "--------------", 
         "--------------", "--------", "----------------");
  for(iter = _Devices.begin(); iter != _Devices.end(); iter++)
  {
    printf("%6d | %14.4f | %14.4f | %8.4f | %16.4f\n",
           (*iter)->GetDeviceID(), (*iter)->GetAverageCPUInnerTime(),
           (*iter)->GetAverageCPUOuterTime(), (*iter)->GetAverageGPUTime(),
           1.0 / (*iter)->GetAverageCPUOuterTime());
  }

  printf("    r%d | %14.4f | %14.4f | %8.4f | %16.4f\n",
         _Renderer->GetDeviceID(), _Renderer->GetAverageCPUInnerTime(),
         _Renderer->GetAverageCPUOuterTime(), _Renderer->GetAverageGPUTime(),
         1.0 / _Renderer->GetAverageCPUOuterTime());
}

void Simulator::Start()
{
  _PauseFlag = false;
  _StopFlag = false;

  // spread parts evenly over all four
  // pthread_barrier_init(&_Barrier, 0, 4);
  // int nPer = _Parameters.N / 4;
  // int extra = _Parameters.N - (nPer*4);
  // _Devices.push_back(new DeviceController(0, 0, nPer-1, &_Barrier));
  // _Devices.push_back(new DeviceController(2, nPer, 2*nPer-1, &_Barrier));
  // _Devices.push_back(new DeviceController(3, 2*nPer, 3*nPer-1, &_Barrier));
  // _Devices.push_back(new DeviceController(4, 3*nPer, 4*nPer-1+extra, &_Barrier));

  // spread parts evenly over first two, which have most memory bandwidth
  // pthread_barrier_init(&_Barrier, 0, 2);
  // int nPer = _Parameters.N / 2;
  // int extra = _Parameters.N % 2;
  // _Devices.push_back(new DeviceController(0, 0, nPer-1, &_Barrier));
  // _Devices.push_back(new DeviceController(2, nPer, 2*nPer-1+extra, &_Barrier));

  // just use a single device
  // pthread_barrier_init(&_Barrier, 0, 1);
  // _Devices.push_back(new DeviceController(0, 0, _Parameters.N - 1, &_Barrier));

  // Three devices for simulation and one for rendering:
  pthread_barrier_init(&_Barrier, 0, 4);
  int nPer = _Parameters.N / 3;
  int extra = _Parameters.N - (nPer*3);
  _Devices.push_back(new DeviceController(2, 0, nPer-1+extra, &_Barrier));
  _Devices.push_back(new DeviceController(3, nPer+extra, 2*nPer-1+extra, &_Barrier));
  _Devices.push_back(new DeviceController(4, 2*nPer+extra, 3*nPer-1+extra, &_Barrier));
  _Renderer = new Renderer(0, &_Barrier, _Vis);
}

void Simulator::Stop()
{
  _PauseFlag = false;
  _StopFlag = true;

  // TODO: wait for devices to stop, then cleanup
  // std::vector<DeviceController*>::iterator iter;
  // for(iter = _Devices.begin(); iter != _Devices.end(); iter++)
  //   delete (*iter);
  // _Devices.clear();
  // TODO: cleanup barrier
  // TODO: cleanup renderer

}

void Simulator::Pause()
{
  _PauseFlag = true;
}

void Simulator::Resume()
{
  _PauseFlag = false;
}

void Simulator::PrintParticles()
{
  for(int i = 0; i < _Parameters.N; i++)
  {
    string format = "[%02i]:  ";
    format += "X=%10.4g, Y=%10.4g,  ";
    format += "VX=%10.4g, VY=%10.4g  ";
    format += "AX=%10.4g, AY=%10.4g\n";
    printf(format.c_str(), i, 
	   _HostX[i], _HostY[i], 
	   _HostVX[i], _HostVY[i], 
	   _HostAX[i], _HostAY[i]);
  }
}
